/*------------------------------------------------------------
This file runs the analyses underlying the tables and figures
in the paper and provides the needed output for use in the paper.  
------------------------------------------------------------*/
capture log close
set more off
set matsize 10000
set scheme s1mono
clear

	*** Set path of working directory here:
		di "Enter working directory: " _request(path) 
		local path $path
	
	cd "`path'"

* Read in data
use "Data/data_for_main_analyses.dta", clear


* Variable labels
* Matters especially for descriptive stats table
label variable saleprice_2010USD "Sale price (2010USD)"
label variable flood_07 "Old floodplain"
label variable post_BW12 "Post Biggert-Waters"
label variable flood07XpostBW12 "Old floodplain*Post Biggert-Waters"
label variable flood_sandy "Flooded by Sandy"
label variable post_Sandy "Post Sandy"
label variable floodSandyXpostSandy "Flooded by Sandy*Post Sandy"
label variable alltreated "New floodplain"
label variable post_ABFE "Post ABFE"
label variable post_prelimWorkMaps "Post prelim. work maps"
label variable post_prelimFIRMs "Post prelim. FIRMs"
label variable trtXpost_newmap "New floodplain*post new maps"

* Regression locals
local post "1.post_ABFE 1.post_prelimWorkMaps 1.post_prelimFIRMs"
local sandy "floodsandyXpostsandyX1flood07 floodsandyXpostsandyX0flood07"
local sandy2 "floodsandyXpostsandyX1flood07 floodsandyXpostsandyX0flood07 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth"
local sandy2labels "floodsandyXpostsandyX0flood07 "Sandy*not in old FP" floodsandyXpostsandyX1flood07 "Sandy*in old FP" 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth "Sandy*depth*in old FP" 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth "Sandy*depth*not in old FP" "
local newmaps "trtXpost_newmapX1flood_sandy trtXpost_newmapX0flood_sandy"
local newmaps2 "trtXp_nmapX*fld_sandyX*fld96"
local timecontrols "yrwk" 
local xsec "i.flood_07##i.flood_sandy##i.alltreated"
local xsec2 "i.flood_07##i.flood_sandy##i.flood_13_1"
local hfiaa "i.postHFIAA##i.flood_13_6"	
local hed_cntrls "c.age c.landsquarefeet##c.landsquarefeet##c.landsquarefeet c.grosssquarefeet##c.grosssquarefeet##c.grosssquarefeet  i.res_units i.com_units i.tot_units i.bld_cls"


		
*** Main Paper***
* Figures
*** Figure 1 - Output from ArcGIS

*** Figure 2 - BW raw version
	* Locals for vertical lines
	local BW = td(10jul2012)
	local Sandy = td(30oct2012)  
	local ABFE = td(28jan2013)
	local workmaps = td(10jun2013)
	local prelimfirms = td(30jan2015)
	local deg "0"
quietly reghdfe logsaleprice_2010USD  ///
	if nbrhdFEsample==1 & blockFEsample==1, absorb(i.block ) residuals(resid)
twoway lpoly resid saledate if flood_07==0, degree(`deg') || ///
	lpoly resid saledate if flood_07==1, degree(`deg') lpattern(dash) || ///
	, xtitle("") ytitle("Residual log sale price"  "", margin(b+5)) ///
	xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
	text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
	xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
	xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
	legend(label(1 "Not in old floodplain") label(2 "In old floodplain") rows(1) region(lwidth(none)))
graph export "output_materials/Graphs/DD_BW_raw.eps", as(eps) replace


*** Figure 3 - Sandy raw version
twoway lpoly resid saledate if flood_sandy==0, degree(`deg') || ///
	lpoly resid saledate if flood_07==1 & flood_sandy==1, degree(`deg') lpattern(dash) || ///
	lpoly resid saledate if flood_07==0 & flood_sandy==1, degree(`deg') lpattern(dot) lcolor(black) || ///
	, xtitle("") ytitle("Residual log sale price" "", margin(b+5)) ///
	xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
	text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `ABFE' "ABFE", place(w) orientation(vertical) size(vsmall)) ///
	xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
	xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
	legend(label(1 "Not flooded") label(2 "Flooded, in old FP") label(3 "Flooded, not in old FP") rows(1) region(lwidth(none)))
graph export "output_materials/Graphs/DD_Sandy_raw.eps", as(eps) replace

*** Figure 4 - Maps raw version

twoway lpoly resid saledate if alltreated==0, degree(`deg') || ///
	lpoly resid saledate if flood_sandy==1 & alltreated==1, degree(`deg') lpattern(dash) || ///
	lpoly resid saledate if flood_sandy==0 & alltreated==1, degree(`deg') lpattern(dot) lcolor(black) || ///
	, xtitle("") ytitle("Residual log sale price" "") ///
	xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) ///
	xline(`ABFE', lpattern(dot) lcolor(gs6)) xline(`workmaps', lpattern(dot) lcolor(gs6)) xline(`prelimfirms', lpattern(dot) lcolor(gs6)) ///
	text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `workmaps' "Prelim. work maps", place(w) orientation(vertical) size(vsmall)) ///
	text(0 `prelimfirms' "Prelim. FIRMs", place(w) orientation(vertical) size(vsmall)) ///
	xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
	xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
	legend(label(1 "Not in new FP") label(2 "In new FP, flooded") label(3 "In new FP, not flooded") rows(1) region(lwidth(none)))
graph export "output_materials/Graphs/DD_FPmaps_raw.eps", as(eps) replace

	

*** Figure 5 - Heterogeneity by Structure Value
* Calculating avg obs per cluster 
	preserve
	keep if nbrhdFEsample==1 & blockFEsample==1
	collapse (count) obs_within_cluster=logsaleprice_2010USD, by(census_tract_INT)
	display "Obs per cluster, columns 1-2"
	summarize obs_within_cluster
	restore
	preserve
	keep if lotFEsample==1
	collapse (count) obs_within_cluster=logsaleprice_2010USD, by(census_tract_INT)
	display "Obs per cluster, columns 3-4"
	summarize obs_within_cluster
	restore


* merge assessment in year of sale onto sale observation by BBLE and year
	gen BBLE=string(block*10000+lot+borough*1000000000,"%10.0f")
	***Use only 2012 assesment data with all sales data
		capture drop _m
		rename saleyear year
		merge m:1 BBLE   using "Data\intermediate_data_files\TC1_Assessments_2012_MKT_VALUES.dta"
		rename year saleyear
		
* Generate Binned Variables Based on structure values		
		drop str_value // based on assessment "market value"
		gen structure_value= (1-lnd_val_ratio)*saleprice

		capture drop str_val_100k_bin BW12_SV_treated sandy_SV_treated maps_SV_treated
		
		egen str_val_100k_bin=cut(structure_value), at(0(100000)500000)
		replace str_val_100k_bin=500000 if str_val_100k_bin==.


		replace str_val_100k_bin=1 if str_val_100k_bin==0

		gen BW12_SV_treated= flood_07*str_val_100k_bin
		gen sandy_SV_treated= flood_sandy*str_val_100k_bin
		gen maps_SV_treated= trtXpost_newmap*str_val_100k_bin


* Mo, yr, Block FEs, big sample
	reghdfe saleprice_2010USD i.flood_07##i.flood_sandy##i.alltreated##ibn.str_val_100k_bin i.flood_07 i.flood_07#ibn.str_val_100k_bin i.post_BW12#ibn.str_val_100k_bin i.post_BW1#ibn.BW12_SV_treated ///
									i.flood_sandy i.flood_sandy#ibn.str_val_100k_bin i.post_Sandy#ibn.str_val_100k_bin  i.post_Sandy#1.flood_07#ibn.sandy_SV_treated i.post_Sandy#0.flood_07#ibn.sandy_SV_treated 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth#ibn.str_val_100k_bin 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth#ibn.str_val_100k_bin  ///
									alltreated i.alltreated#ibn.str_val_100k  i.post_ABFE#ibn.str_val_100k ibn.flood_sandy#ib0.maps_SV_treated     ///
									if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT str_val_100k_bin yrwk) vce(cluster census_tract_INT)
				
		estimates store A

						coefplot  (A, keep(1.flood_sandy#*maps_SV_treated) label("In New Floodplain - Sandy Flooded")) (A, keep(0.flood_sandy#*maps_SV_treated)  label("In New Floodplain - Not Flooded") ) , bylabel("Block_FE") || ///
								,vert  yline(0, lwidth(vthin)) xline(3.5, lcolor(gray) lwidth(vthin) lpattern(dash)) xline(  6.5, lwidth(medthick)) xline(9.5, lcolor(gray) lwidth(vthin) lpattern(dash)) /// 
								xtitle("Estimated Structure Value - $100k Bins")   ///
								coeflabel(	1.flood_sandy#1.maps_SV_treated = "$0-$100k"  		///
								1.flood_sandy#100000.maps_SV_treated = "$100-$200k"  	///
								1.flood_sandy#200000.maps_SV_treated = "$200-$300k"  	///
								1.flood_sandy#300000.maps_SV_treated = "$300-$400k"  	///
								1.flood_sandy#400000.maps_SV_treated = "$400-$500k"  	///
								1.flood_sandy#500000.maps_SV_treated = "Over $500"  	///
							0.flood_sandy#1.maps_SV_treated = "$0-$100k"  		///
								0.flood_sandy#100000.maps_SV_treated = "$100-$200k"  	///
								0.flood_sandy#200000.maps_SV_treated = "$200-$300k"  	///
								0.flood_sandy#300000.maps_SV_treated = "$300-$400k"  	///
								0.flood_sandy#400000.maps_SV_treated = "$400-$500k"  	///
								0.flood_sandy#500000.maps_SV_treated = "Over $500" ) 	///
								xlabel(,angle(35) ) ///
								legend(rows(1) pos(12) bmargin(off) region(lstyle(none)))  byopts(row(1) title("Coef Estimates by Structure Value"))	 ///
								ytitle("Implied Change in Property Value") 
				graph export "output_materials/Graphs/DD_FPmaps_Hetero_by_Str_Val.eps", as(eps) replace

	
* Tables
* Read in data
use "Data/data_for_main_analyses.dta", clear	

*** Table 1
		eststo clear
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'  if nbrhdFEsample==1 & blockFEsample==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'  if lotFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Lot FEs -  
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps'  if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
		* Table
		esttab using "output_materials/Tables/DDmodels_heterog_dmgcontrol.tex", keep(1.flood_07#1.post_BW12 *floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy   ) /// 
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			mtitles("Neighborhood FE" "Block FE" "Block FE" "Lot FE") nonotes label ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace		

			
		
			

*** Table 2 - Alternative Specifications ***
			eststo clear

			*	1.) Day-of-Sale FE (robustness to flexible time trends)
						local newmaps "trtXpost_newmapX1flood_sandy trtXpost_newmapX0flood_sandy"
						local xsec "i.flood_07##i.flood_sandy##i.alltreated"
						local sandy2 "floodsandyXpostsandyX1flood07 floodsandyXpostsandyX0flood07 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth"

					local timecontrols "saledate" 
					
					*** Block FEs -  
					eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 i.flood_sandy `sandy2' alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
					gen date_fe_sample = e(sample)
					*** Lot FEs -  
					eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
			
			
			
			
			*	2.) Borough*yr-mo FE (robustness to flexible time trends)
					local newmaps "trtXpost_newmapX1flood_sandy trtXpost_newmapX0flood_sandy"
					local xsec "i.flood_07##i.flood_sandy##i.alltreated"
					local sandy2 "floodsandyXpostsandyX1flood07 floodsandyXpostsandyX0flood07 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth"

					local timecontrols  "boro_yrmo"
					
					*** Block FEs -  
					eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 i.flood_sandy `sandy2' alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
					gen boro_mo_fe_sample = e(sample)
					*** Lot FEs -  
					eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)

						
				
				
			*	3.) Quadratic in Sandy inundation 
					* Regression locals
					local newmaps "trtXpost_newmapX1flood_sandy trtXpost_newmapX0flood_sandy"
					local timecontrols "yrwk" 
					local xsec "i.flood_07##i.flood_sandy##i.alltreated"
					
					* Redefine Sandy2 variable to include the squared terms in depth
					local sandy2 "floodsandyXpostsandyX1flood07 floodsandyXpostsandyX0flood07 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth#c.sandy_inund_depth 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth#c.sandy_inund_depth"
					local sandy2labels "floodsandyXpostsandyX0flood07 "Sandy*not in old FP" floodsandyXpostsandyX1flood07 "Sandy*in old FP" 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth "Sandy*depth*in old FP" 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth "Sandy*depth*not in old FP" 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth#c.sandy_inund_depth "Sandy*depth^2*in old FP" 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth#c.sandy_inund_depth "Sandy*depth^2*not in old FP" "
				
					*** Block FEs -  
					eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 i.flood_sandy `sandy2' alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
					*** Block FEs -  

					*** Lot FEs -  
					eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)


			* Output the Table
					local sandy2labels "floodsandyXpostsandyX0flood07 "Sandy*not in old FP" floodsandyXpostsandyX1flood07 "Sandy*in old FP" 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth "Sandy*depth*in old FP" 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth "Sandy*depth*not in old FP" 1.floodsandyXpostsandyX1flood07#c.sandy_inund_depth#c.sandy_inund_depth "Sandy*depth\$^2$*in old FP" 1.floodsandyXpostsandyX0flood07#c.sandy_inund_depth#c.sandy_inund_depth "Sandy*depth\$^2$*not in old FP" "

				esttab using "output_materials/Tables/DDmodels_Alt_Specs_Table.tex", 	keep(1.flood_07#1.post_BW12 floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth* trtXpost_newmapX*flood_sandy) ///
					label numbers nonotes  collabels(none)  ///
					coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
					mgroups("\underline{Sale Date FE}"  "\underline{Boro*Yr-Mo FE}"  "\underline{Depth-Squared}", pattern(1 0 1 0 1 0) span prefix(\multicolumn{2}{c}{) suffix(})) ///
					mtitles( "Block FE" "Lot FE" "Block FE" "Lot FE" "Block FE" "Lot FE" ) ///
					se star(* 0.10 ** 0.05 *** 0.01) compress tex replace
			
	
	
	

*** Table 3 - Matched D-in-Ds
		preserve
					
				gen nm_post=trtXpost_newmap
				replace nm_post=1 if post_ABFE==1 & flood_13_1==0 & flood_13_6==0 & flood_15==0
				
						

		*** Calculate Neighborhood, Block, and Lot-specific pre-treatment mean price
				keep if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1

						
				gen pre_price=saleprice_2010USD if post_BW12==0

				egen nbrhd_pre_price=mean(pre_price), by(neighborhoodINT)		
				egen block_pre_price=mean(pre_price), by(blockINT)		
				egen lot_pre_price=mean(pre_price), by(lotINT)

				foreach geo in nbrhd block lot {
				gen `geo'_ln_pre_price=ln(`geo'_pre_price)

				}
				

		*** Calculate Neighborhood, Block, and Lot-specific pre-treatment mean price trends
				
		sort lotINT saledate
		* Crude previous sale trend
			* Calculate implied mean annual price change between each sale in the pre-period
				bysort lotINT: generate aagr = ((pre_price - pre_price[_n-1]) / (pre_price[_n-1]*(saledate - saledate[_n-1])))*365
				
			* Take average of pre-treatment lot-specific price change rates over N, B, & L
				bysort neighborhoodINT: egen nbrhd_avgaagr_preBW=mean(aagr)
				bysort blockINT: egen block_avgaagr_preBW=mean(aagr)
				bysort lotINT: egen lot_avgaagr_preBW = mean(aagr) 


		*** Generate Indicator Variables for Fixed Effects (diff command doesn't accept factor variables)

			*** Week of sample FE
				sum yrwk
				local max_yrwk=r(max)
				qui levelsof yrwk, local(levels)
				
				display "This will take a moment"
				
				foreach l of local levels {
					qui gen fe_`l'= (yrwk==`l')
					*display "`l' of `max_yrwk'"
				}	

		


		eststo clear
		*** do D-in-D with matching
		eststo: diff logsaleprice_2010USD if nbrhdFEsample==1 & blockFEsample==1 & flood_sandy==0, period(nm_post)  treated(alltreated) kernel id(neighborhoodINT) cov( grosssquarefeet landsquarefeet age nbrhd_pre_price lat lon block_avgaagr_preBW ) addcov(grosssquarefeet landsquarefeet age  nbrhd_pre_price block_avgaagr_preBW fe_* )  cluster(census_tract_INT ) robust support
				 
				 
		eststo: diff logsaleprice_2010USD if nbrhdFEsample==1 & blockFEsample==1  & flood_sandy==0, period(nm_post) treated(alltreated) kernel id(blockINT) cov( grosssquarefeet landsquarefeet age block_pre_price lat lon block_avgaagr_preBW ) addcov(grosssquarefeet landsquarefeet age   block_pre_price block_avgaagr_preBW fe_* )  cluster(census_tract_INT ) robust  support
				
		eststo: diff logsaleprice_2010USD if lotFEsample==1 & blockFEsample==1  & flood_sandy==0, period(nm_post) treated(alltreated) kernel id(blockINT) cov( grosssquarefeet landsquarefeet age block_pre_price lat lon block_avgaagr_preBW ) addcov(grosssquarefeet landsquarefeet age block_pre_price  nbrhd_pre_price block_avgaagr_preBW fe_* )  cluster(census_tract_INT ) robust support
				
		eststo: diff logsaleprice_2010USD if blockFEsample==1 & lotFEsample==1  & flood_sandy==0, period(nm_post) treated(alltreated) kernel id(lotINT) cov( grosssquarefeet landsquarefeet age lot_pre_price lat lon block_avgaagr_preBW ) addcov(grosssquarefeet landsquarefeet age  lot_pre_price block_avgaagr_preBW fe_*) cluster(census_tract_INT ) robust support

		 label var _diff "Floodplain maps"
		
			esttab using "output_materials/Tables/matched_did.tex", 	keep(_diff) ///
 			label numbers nonotes  collabels(none)  ///
			coeflabels(_diff "Floodplain maps") ///
			mgroups( "\underline{Pre-Treatment Price Spatial Aggregation by:}"  , pattern(1 0 0 0) span prefix(\multicolumn{4}{c}{) suffix(})) ///
			mtitles( "Neighborhood" "Block" "Block" "Lot") ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace		
restore
			

*** Table 4 - Alternative Samples ***

			* 1.) Exclude buffer both inside and outside of FP from sample to address possible spillovers
					gen exclude= dist_2007_FP_edge<50
					eststo clear


					*** Block FEs -  
					eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 i.flood_sandy `sandy2' alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1 & exclude!=1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
					*** Lot FEs -  
					eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1 & exclude!=1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)

					
				
			* 2.)	Use donut outside floodplain as control group 
					capture drop exclude
					gen exclude= dist_2007_FP_edge>500 & flood_07==0
				

					*** Block FEs -  
					eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 i.flood_sandy `sandy2' alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1 & exclude!=1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
					*** Lot FEs -  
					eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1 & exclude!=1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
			
			
			* 3.) Rerun main estimates on single-family homes
				*** Block FEs -  
				eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1 & residentialunits==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
				*** Lot FEs -  
				eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps'  if lotFEsample==1 & residentialunits==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
			
			
			
			
				
			*Output the Table
				esttab using "output_materials/Tables/DDmodels_Alt_Samples_Table.tex", 	keep(1.flood_07#1.post_BW12 floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth* trtXpost_newmapX*flood_sandy) ///
					label numbers nonotes  collabels(none)  ///
					coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
					mgroups( "\underline{Exclude 50M Boundary}"     "\underline{500M Buffer as Cntrl}"  "\underline{Single Family Homes}" , pattern(1 0 1 0 1 0) span prefix(\multicolumn{2}{c}{) suffix(})) ///
					mtitles( "Block FE" "Lot FE" "Block FE" "Lot FE" "Block FE" "Lot FE" "Block FE" "Lot FE") ///
					se star(* 0.10 ** 0.05 *** 0.01) compress tex replace			
			
	

	

*** Table 5 - Census Tract Comparison 
		preserve
			use "Data\ACS_11_16_Diffs.dta", clear
			gen obs_count=1
			
			keep in_fp diff_11_16_* se_diff_* obs_count
			
			qui ds
			local var_order=r(varlist)
			 
			
			collapse (mean) diff_11_16_*  (sum) se_diff_* obs_count, by(in_fp)
			
			order `var_order'
			foreach var of varlist se_* {
				replace `var'= sqrt(`var'/obs_count^2)
			}


			rename diff_11_16_* diff_*
			xpose, clear varname
			gen se= (substr(_varname,1,3)=="se_")
			replace _varname=substr(_varname,4,length(_varname)-3) if substr(_varname,1,3)=="se_"
			egen var_numb=group(_varname)
			
			rename v1 out_of_fp_mean
			rename v2 in_fp_mean

			drop if _varname=="in_fp"
			
			reshape wide out_of_fp_mean in_fp_mean, i(var_numb) j(se)
			
			rename out_of_fp_mean0 out_of_fp_mean
			rename out_of_fp_mean1 out_of_fp_se 
			
			rename in_fp_mean0 in_fp_mean
			rename in_fp_mean1 in_fp_se 
			
			gen diff= out_of_fp_mean-in_fp_mean
			gen diff_se= sqrt(out_of_fp_se^2 + in_fp_se^2)
			gen diff_tcalc= diff/diff_se
			
			gen diff_sig=string(round(diff_tcalc,.001))
			replace diff_sig=diff_sig+"*" if abs(diff_tcalc)>1.645
			replace diff_sig=diff_sig+"*" if abs(diff_tcalc)>1.960
			replace diff_sig=diff_sig+"*" if abs(diff_tcalc)>2.575

		drop var_num

		order _varname out_of_fp_mean in_fp_mean out_of_fp_se in_fp_se _varname diff diff_se diff_tcalc diff_sig

		export excel using"output_materials\Tables\ACS_Balance_Table.xlsx", replace firstrow(variables)

		restore



*** Table 6 - Tests for Sample Selection	
		* Uses pre-treatment sale prices to test for selection into post-treatment transactions
		* Variable handling
		gen pretreat = (saledate < td(01jun2012))
		local land "c.landsquarefeet##c.landsquarefeet##c.landsquarefeet##c.landsquarefeet"
		local floor "c.grosssquarefeet##c.grosssquarefeet##c.grosssquarefeet##c.grosssquarefeet"
		local yearbuilt "c.yearbuilt##c.yearbuilt##c.yearbuilt##c.yearbuilt"
		local units "c.residentialunits##c.residentialunits##c.residentialunits##c.residentialunits"
		reghdfe logsaleprice_2010USD `land' `floor' `yearbuilt' `units' if nbrhdFEsample==1 & blockFEsample==1 & pretreat==1, absorb(neighborhoodINT, savefe) resid
		predict selectionindex if nbrhdFEsample==1 & blockFEsample==1


		* Regressions
		* Note lot FE impossible because avg. pre-treatment sale price doesn't vary within lot
		eststo clear
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe selectionindex `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' if nbrhdFEsample==1 & blockFEsample==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe selectionindex `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe selectionindex `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Lot FEs -  
		eststo: reghdfe selectionindex i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
		* Table
		esttab using "output_materials/Tables/Selection.tex", keep(1.flood_07#1.post_BW12 floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy) ///
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			mtitles("Neighborhood FE" "Block FE" "Block FE" "Lot FE") nonotes ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace


	
	
	
*** Appendix***	

* Read in data
use "Data/data_for_main_analyses.dta", clear

* Figures	

*** Figure A1: Sales Counts by Year and Borough
	set scheme s1mono
		preserve
			gen count=1

			collapse (sum) count, by(saleyear borough)
				drop if borough==.
			reshape wide count, i(saleyear) j(borough)

			rename count1 Manhattan
			rename count2 Bronx
			rename count3 Brooklyn
			rename count4 Queens
			rename count5 Staten_Island
			
			

			graph bar (sum) Manhattan Bronx Brooklyn Queens Staten_Island , over(saleyear, label(angle(45))) stack legend( label(1 Manhattan) label(2 Bronx) ///
						label(3 Brooklyn) label(4 Queens) label(5 "Staten Island") angle(45)) ytitle("Number of Transactions by Year") graphregion(color(white))

			graph export "output_materials/Graphs/sales_counts_by_year.eps", as(eps) replace
			restore	
		
		


***	Figure A2: Venn Diagram of Treatment Groups

		*** Neighborhood FE Sample	
			quiet{ 
				count if flood_07==1 & nbrhdFEsample==1 & blockFEsample==1
			local old_FP=r(N)
			
			count if flood_sandy==1 & nbrhdFEsample==1 & blockFEsample==1
			local sandy_flooded=r(N)
			
			count if alltreated==1 & nbrhdFEsample==1 & blockFEsample==1
			local new_FP=r(N)
			
			
			count if flood_07==1 & flood_sandy==1 & nbrhdFEsample==1 & blockFEsample==1
			local old_FP_sandy=r(N)
			
			
			count if flood_07==1 & alltreated==1 & nbrhdFEsample==1 & blockFEsample==1
			local old_FP_new_FP=r(N)
			
			count if flood_sandy==1 & alltreated==1 & nbrhdFEsample==1 & blockFEsample==1
			local sandy_new_FP=r(N)
			
			count if flood_07==1 & flood_sandy==1 & alltreated==1 & nbrhdFEsample==1 & blockFEsample==1
			local in_all_3=r(N)
			
			
			count if nbrhdFEsample==1 & blockFEsample==1
			local in_sample=r(N)
				pvenn2, pA(`old_FP') pB(`sandy_flooded') pAB(`old_FP_sandy') 	pC(`new_FP') pAC(`old_FP_new_FP') pBC(`sandy_new_FP')  ang(30)  ptotal(`in_sample') drawlabel(1) lc(gs9 gs4 gs1 gs1) lp(solid dash dot solid )  lw(thin thin thin medium) yscale(lstyle(none)) xscale(lstyle(none)) plotregion(color(white) lwidth(none)) legend(order(1 "Old FP" 3 "Sandy Flooded" 5 "New FP" 7 "Nbrhd FE Sample")) 
			}

		*** Repeated sales sample - Figure not presented
			quiet{ 
						count if flood_07==1 & lotFEsample==1
					local old_FP=r(N)
					
					count if flood_sandy==1 & lotFEsample==1
					local sandy_flooded=r(N)
					
					count if alltreated==1 & lotFEsample==1
					local new_FP=r(N)
					
					
					count if flood_07==1 & flood_sandy==1 & lotFEsample==1
					local old_FP_sandy=r(N)
					
					
					count if flood_07==1 & alltreated==1 & lotFEsample==1
					local old_FP_new_FP=r(N)
					
					count if flood_sandy==1 & alltreated==1 & lotFEsample==1
					local sandy_new_FP=r(N)
					
					count if flood_07==1 & flood_sandy==1 & alltreated==1 & lotFEsample==1
					local in_all_3=r(N)
					
					
					count if lotFEsample==1
					local in_sample=r(N)
						pvenn2, pA(`old_FP') pB(`sandy_flooded') pAB(`old_FP_sandy') 	pC(`new_FP') pAC(`old_FP_new_FP') pBC(`sandy_new_FP')  pABC(`in_all_3') ang(30)  ptotal(`in_sample') drawlabel(1) lc(gs9 gs4 gs1 gs1) lp(solid dash dot solid )  lw(thin thin thin medium) yscale(lstyle(none)) xscale(lstyle(none)) plotregion(color(white) lwidth(none)) legend(order(1 "Old FP" 3 "Sandy Flooded" 5 "New FP" 7 "Repeated Sales Sample"))
					}
				
*** Figures A3 - A9: Maps Output by ArcGIS

*** Figure A10 - BW raw - diffs
		local BW = td(10jul2012)
		local Sandy = td(30oct2012)  
		local ABFE = td(28jan2013)
		local workmaps = td(10jun2013)
		local prelimfirms = td(30jan2015)
		local deg "0"
		capture drop resid
		quietly reghdfe logsaleprice_2010USD  ///
	if nbrhdFEsample==1 & blockFEsample==1, absorb(i.block ) residuals(resid)

		preserve
		collapse (mean) resid, by(saledate flood_07)
		reshape wide resid, i(saledate) j(flood_07)
		lookfor resid
		generate rdiff = resid1 - resid0
		twoway lpolyci rdiff saledate, degree(`deg') ciplot(rline) || ///
			, xtitle("") ytitle("Diff. in resid. log sale price"  "", margin(b+5)) ///
			xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
			yline(0, lpattern(solid) lcolor(gs15)) ///
			text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
			xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
			xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
			legend(off)
		graph export "output_materials/Graphs/DD_BW_raw_diff.eps", as(eps) replace	
		restore


*** Figure A11 - Sandy raw - diffs
		preserve
		generate group = .
		replace group = 0 if flood_sandy==0
		replace group = 1 if flood_07==1 & flood_sandy==1
		replace group = 2 if flood_07==0 & flood_sandy==1
		drop if missing(group)
		collapse (mean) resid, by(saledate group)
		reshape wide resid, i(saledate) j(group)
		generate rdiff1 = resid1 - resid0
		generate rdiff2 = resid2 - resid0
		twoway lpolyci rdiff1 saledate, degree(`deg') ciplot(rline) || ///
			lpolyci rdiff2 saledate, degree(`deg') lpattern(dash) ciplot(rline) || ///
			, xtitle("") ytitle("Diff. in resid. log sale price" "", margin(b+5)) ///
			xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
			yline(0, lpattern(solid) lcolor(gs15)) ///
			text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `ABFE' "ABFE", place(w) orientation(vertical) size(vsmall)) ///
			xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
			xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
			legend(label(2 "Flooded, in old FP - control") label(4 "Flooded, not in old FP - control") region(lwidth(none)))
		graph export "output_materials/Graphs/DD_Sandy_raw_diff.eps", as(eps) replace
		restore

*** Figure A12 - Maps raw - diffs
		preserve
		generate group = .
		replace group = 0 if alltreated==0
		replace group = 1 if flood_sandy==1 & alltreated==1
		replace group = 2 if flood_sandy==0 & alltreated==1
		drop if missing(group)
		collapse (mean) resid, by(saledate group)
		reshape wide resid, i(saledate) j(group)
		generate rdiff1 = resid1 - resid0
		generate rdiff2 = resid2 - resid0
		twoway lpolyci rdiff1 saledate, degree(`deg') ciplot(rline) || ///
			lpolyci rdiff2 saledate, degree(`deg') lpattern(dash) ciplot(rline) || ///
			, xtitle("") ytitle("Diff. in resid. log sale price" "") ///
			xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) ///
			xline(`ABFE', lpattern(dot) lcolor(gs6)) xline(`workmaps', lpattern(dot) lcolor(gs6)) xline(`prelimfirms', lpattern(dot) lcolor(gs6)) ///
			yline(0, lpattern(solid) lcolor(gs15)) ///
			text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `workmaps' "Prelim. work maps", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `prelimfirms' "Prelim. FIRMs", place(w) orientation(vertical) size(vsmall)) ///
			xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
			xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
			legend(label(2 "In new FP, flooded - control") label(4 "In new FP, not flooded - control") region(lwidth(none)))
		graph export "output_materials/Graphs/DD_FPmaps_raw_diff.eps", as(eps) replace
		restore	
		


*** Figure A13 - BW raw version
	
		* Locals for vertical lines
		local BW = td(10jul2012)
		local Sandy = td(30oct2012)  
		local ABFE = td(28jan2013)
		local workmaps = td(10jun2013)
		local prelimfirms = td(30jan2015)

		* Lpoly parameters
		local deg "0"
		local bw "120"

		capture drop resid
		quietly reghdfe logsaleprice_2010USD  ///
			if lotFEsample==1, absorb(i.lotINT) residuals(resid)
		twoway lpoly resid saledate if flood_07==0, degree(`deg') bwidth(`bw') || ///
			lpoly resid saledate if flood_07==1, degree(`deg') bwidth(`bw') lpattern(dash) || ///
			, xtitle("") ytitle("Residual log sale price"  "", margin(b+5)) ///
			xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
			text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
			xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
			xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
			legend(label(1 "Not in old floodplain") label(2 "In old floodplain") rows(1) region(lwidth(none)))
		graph export "output_materials/Graphs/DD_BW_raw_w_lot_fe.eps", as(eps) replace


*** Figure A14 - Sandy raw version
		twoway lpoly resid saledate if flood_sandy==0, degree(`deg') bwidth(`bw') || ///
			lpoly resid saledate if flood_07==1 & flood_sandy==1, degree(`deg') bwidth(`bw') lpattern(dash) || ///
			lpoly resid saledate if flood_07==0 & flood_sandy==1, degree(`deg') bwidth(`bw') lpattern(dot) lcolor(black) || ///
			, xtitle("") ytitle("Residual log sale price" "", margin(b+5)) ///
			xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
			text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `ABFE' "ABFE", place(w) orientation(vertical) size(vsmall)) ///
			xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
			xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
			legend(label(1 "Not flooded") label(2 "Flooded, in old FP") label(3 "Flooded, not in old FP") rows(1) region(lwidth(none)))
		graph export "output_materials/Graphs/DD_Sandy_raw_w_lot_fe.eps", as(eps) replace



*** Figure A15 - Maps raw version
		twoway lpoly resid saledate if alltreated==0, degree(`deg') bwidth(`bw') || ///
			lpoly resid saledate if flood_sandy==1 & alltreated==1, degree(`deg') bwidth(`bw') lpattern(dash) || ///
			lpoly resid saledate if flood_sandy==0 & alltreated==1, degree(`deg') bwidth(`bw') lpattern(dot) lcolor(black) || ///
			, xtitle("") ytitle("Residual log sale price" "") ///
			xline(`BW', lpattern(dot) lcolor(gs6)) xline(`Sandy', lpattern(dot) lcolor(gs6)) ///
			xline(`ABFE', lpattern(dot) lcolor(gs6)) xline(`workmaps', lpattern(dot) lcolor(gs6)) xline(`prelimfirms', lpattern(dot) lcolor(gs6)) ///
			text(0 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `workmaps' "Prelim. work maps", place(w) orientation(vertical) size(vsmall)) ///
			text(0 `prelimfirms' "Prelim. FIRMs", place(w) orientation(vertical) size(vsmall)) ///
			xtick(15706 16071 16437 16802 17167 17532 17898 18263 18628 18993 19359 19724 20089 20454 20820, angle(0)) ///
			xlabel(15706 "2003" 16071 "2004" 16437 "2005" 16802 "2006" 17167 "2007" 17532 "2008" 17898 "2009" 18263 "2010" 18628 "2011" 18993 "2012" 19359 "2013" 19724 "2014" 20089 "2015" 20454 "2016" 20820 "2017", angle(45) notick labgap(2)) ///
			legend(label(1 "Not in new FP") label(2 "In new FP, flooded") label(3 "In new FP, not flooded") rows(1) region(lwidth(none)))
		graph export "output_materials/Graphs/DD_FPmaps_raw_w_lot_fe.eps", as(eps) replace


*** Figure A16 - Google Trends
	preserve
					
			* Reading data
			use "Data/Google trends/Google_trends_master.dta", clear

			* Residualizing on month dummies and graphing
			regress search i.month if area=="New York NY" & term=="floodplain"
			predict r1 if e(sample), residuals

			regress search i.month if area=="United States" & term=="floodplain"
			predict r2 if e(sample), residuals
			local Katrina = 547
			local BW = 630
			local Sandy = 633
			local ABFE = 636
			local workmaps = 641
			local prelimfirms = 660

			* Combined graph
			twoway line r1 modate if area=="New York NY" & term=="floodplain" & inrange(modate, tm(2011m1), tm(2016m12)) ///
				|| line r2 modate if area=="United States" & term=="floodplain" & inrange(modate, tm(2011m1), tm(2016m12)), lpattern(dash) ///
				xline(`Sandy', lpattern(dot) lcolor(gs6)) /* Oct 2012 Sandy */ ///
				xline(`BW', lpattern(dot) lcolor(gs6)) /* July 2012 BW */ ///
				xline(`ABFE', lpattern(dot) lcolor(gs6)) /* Jan 2013 ABFE */ ///
				xline(`workmaps', lpattern(dot) lcolor(gs6)) /* Prelim work maps June 2013 */ ///
				xline(`prelimfirms', lpattern(dot) lcolor(gs6)) /* Prelim. FIRMs Jan 2015 */ ///
				text(-15 `BW' "Biggert Waters Act", place(w) orientation(vertical) size(vsmall)) ///
				text(-15 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
				text(-15 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
				text(-15 `workmaps' "Prelim. work maps", place(w) orientation(vertical) size(vsmall)) ///
				text(-15 `prelimfirms' "Prelim. FIRMs", place(w) orientation(vertical) size(vsmall)) ///
				legend(label(1 "NYC") label(2 "USA") rows(1) region(lwidth(none)) colgap(*2))
			graph export "output_materials/Graphs/Google_floodplain_searches_comb.eps", replace

	restore

*** Figure A17 - Plot NFIP Policies in force in NYC
	preserve		
			* Reading data
			import excel using "Data/NFIP_Policy_Statistics/NYC_NFIP_data.xlsx", sheet("NFIP") firstrow case(lower) clear

			* Variable handling
			format date %td
			generate tm = mofd(date)
			format tm %tm


			***********************
			* Graphs
			***********************

			* Locals for vertical lines
			local Sandy = td(30oct2012) 
			local ABFE = td(28jan2013)
			local workmaps = td(10jun2013)
			local prelimfirms = td(30jan2015)

			twoway scatter policies date ///
				, xline(`Sandy', lpattern(dot) lcolor(gs6)) ///
				xline(`ABFE', lpattern(dot) lcolor(gs6)) ///
				xline(`prelimfirms', lpattern(dot) lcolor(gs6)) ///
				text(45000 `Sandy' "Sandy", place(w) orientation(vertical) size(vsmall)) ///
				text(45000 `ABFE' "ABFE maps", place(w) orientation(vertical) size(vsmall)) ///
				text(45000 `prelimfirms' "Prelim. FIRMs", place(w) orientation(vertical) size(vsmall))
			graph export "output_materials/Graphs/NFIP_policies.eps", as(eps) replace
	restore

*** Figures A18 & A19: Maps Output by ArcGIS


*** Figure A20 - Matched D-in-Ds 

	preserve
			
		gen nm_post=trtXpost_newmap
		replace nm_post=1 if post_ABFE==1 & flood_13_1==0 & flood_13_6==0 & flood_15==0
		
				

	*** Calculate Neighborhood, Block, and Lot-specific pre-treatment mean price
		keep if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1

				
		gen pre_price=saleprice_2010USD if post_BW12==0

		egen nbrhd_pre_price=mean(pre_price), by(neighborhoodINT)		
		egen block_pre_price=mean(pre_price), by(blockINT)		
		egen lot_pre_price=mean(pre_price), by(lotINT)

		foreach geo in nbrhd block lot {
		gen `geo'_ln_pre_price=ln(`geo'_pre_price)

		}
		

	*** Calculate Neighborhood, Block, and Lot-specific pre-treatment mean price trends
			
	sort lotINT saledate
	* Crude previous sale trend
		* Calculate implied mean annual price change between each sale in the pre-period
			bysort lotINT: generate aagr = ((pre_price - pre_price[_n-1]) / (pre_price[_n-1]*(saledate - saledate[_n-1])))*365
		
	* Take average of pre-treatment lot-specific price change rates over N, B, & L
		bysort neighborhoodINT: egen nbrhd_avgaagr_preBW=mean(aagr)
		bysort blockINT: egen block_avgaagr_preBW=mean(aagr)
		bysort lotINT: egen lot_avgaagr_preBW = mean(aagr) 


	*** Generate Indicator Variables for Fixed Effects (diff command doesn't accept factor variables)

	*** Week of sample FE
		sum yrwk
		local max_yrwk=r(max)
		qui levelsof yrwk, local(levels)
		
		display "This will take a moment"
		
		foreach l of local levels {
			qui gen fe_`l'= (yrwk==`l')
			*display "`l' of `max_yrwk'"
		}	

	 

	*** Make Matched D-in-Ds graph	 

		capture drop in_sample				 
		gen in_sample=(saledate < td(01jan2013))		 
						 
			mat A = J(1,14,.)
				mat B = J(1,14,.)
				local i = 0
			 
						 
						 
		forvalues year=2005(1)2012 {
			foreach mon in "jan"  {
					display "`mon'`year'"
					capture drop placebo_treat_date
					gen placebo_treat_date = (saledate < td(01`mon'`year'))
						 
						*** Run placebo treatments at 6 month intervals throughout the pre-period
								diff logsaleprice_2010USD if blockFEsample==1 & lotFEsample==1  & flood_sandy==0 & in_sample==1, period(placebo_treat_date) treated(alltreated) kernel id(lotINT) cov( grosssquarefeet landsquarefeet age lot_pre_price lat lon block_avgaagr_preBW ) addcov(grosssquarefeet landsquarefeet age  lot_pre_price block_avgaagr_preBW fe_*) cluster(census_tract_INT ) robust support
								
				mat A[1,1+`i'] = r(did)
				mat B[1,1+`i'] = r(se_dd)
					local i = `i'+1			
								
							}		 
				}		 
						 
			diff logsaleprice_2010USD if blockFEsample==1 & lotFEsample==1 & flood_sandy==0, period(nm_post) treated(alltreated) kernel id(lotINT) cov( grosssquarefeet landsquarefeet age lot_pre_price lat lon block_avgaagr_preBW ) addcov(grosssquarefeet landsquarefeet age  lot_pre_price block_avgaagr_preBW fe_*) cluster(census_tract_INT ) robust support
			
					mat A[1,1+`i'] = r(did)
				mat B[1,1+`i'] = r(se_dd)
				
					mat colnames A = "jan2007" "jul2007" "jan2008" "jul2008" "jan2009" "jul2009" "jan2010" "jul2010"  "jan2011" "jul2011" "jan2012" "jul2012" "jan2013" "jul2013"
				mat colnames B = "jan2007" "jul2007" "jan2008" "jul2008" "jan2009" "jul2009" "jan2010" "jul2010"  "jan2011" "jul2011" "jan2012" "jul2012" "jan2013" "jul2013"

				estadd mat A = A
				estadd mat B = B
							mat GRAPH_IS = J(9,3,.)
						
							forvalues i=4/9 {
							mat GRAPH_IS[`i'-3,1] = A[1,`i'], A[1,`i'] - 1.96*B[1,`i'], A[1,`i'] + 1.96*B[1,`i']
							}
				 
						 mat rownames  GRAPH_IS =  "Jan2008" "Jan2009" "Jan2010"   "Jan2011" "Jan2012"  "Jan2013" 
						 coefplot mat(GRAPH_IS[,1]), scale(.75)  name(Estimates, replace) subtitle("`1'",color(black)) ///
						graphregion(color(white)) vertical msymbol(O) mcolor(gs6) lwidth(medium) lcolor(gs6) mlwidth(vthin) msize(medium) ///
						xtitle("Assumed Treatment Date") ytitle("Matched D-in-D") yline(0, lcolor(black)) ///
						ci((GRAPH_IS[,2] GRAPH_IS[,3])) connect(direct) ciopts(recast(rcspike) fcolor(gs10) fintensity(50) lcolor(gs10) lwidth(medium))
					graph export "output_materials/Graphs\placebo_mDiD.png", replace
		restore



*** Figure A21: Map Output by ArcGIS



* Appendix Tables *
* Read in data
use "Data/data_for_main_analyses.dta", clear

*** Table A1 - Timeline: Data Entered by Hand

*** Table A2 - Property Counts by Flood Zones

		preserve
				duplicates drop borough block lot easement, force
				count

				replace NFHL_2007 ="A" if NFHL_2007 =="AO" | NFHL_2007 =="AE"

				replace zone_2013_1="X500" if zone_2013_1=="Shaded X"
				replace zone_2013_1="VE" if zone_2013_1=="V"
				replace zone_2013_1="X" if zone_2013_1==""

				replace zone_2013_6="A" if zone_2013_6=="AO" | zone_2013_6=="AE"
				replace zone_2013_6="X500" if zone_2013_6=="Shaded X"
				replace zone_2013_6="X" if zone_2013_6=="Unshaded X"
				replace zone_2013_6="X" if zone_2013_6=="" 


				replace zone_2015="A" if zone_2015=="AO" | zone_2015=="AE"
				replace zone_2015="X500" if zone_2015=="0.2 PCT ANNUAL CHANCE FLOOD HAZARD"
				replace zone_2015="X" if zone_2015==""
				
				tab NFHL_2007 
				tab zone_2013_1
				tab zone_2013_6
				tab zone_2015
		restore



*** Table A3 - Descriptions of flood risk zones

*** Table A4 - Descriptive Stats
		* Descriptive stats - Block FE Sample
		eststo clear
		local variables "saleprice_2010USD flood_07 post_BW12 flood07XpostBW12 flood_sandy post_Sandy floodSandyXpostSandy alltreated post_ABFE post_prelimWorkMaps post_prelimFIRMs trtXpost_newmap"
		estpost summarize `variables' if nbrhdFEsample==1 & blockFEsample==1, listwise
		esttab using "output_materials/Tables/Descriptive_stats.tex", cells("sum mean(fmt(2)) sd min max") substitute(sum Count mean Mean sd Stdev min Min max Max) ///
			label nomtitle nonumber replace

		* Descriptive stats - Lot FE Sample
		eststo clear
		local variables "saleprice_2010USD flood_07 post_BW12 flood07XpostBW12 flood_sandy post_Sandy floodSandyXpostSandy alltreated post_ABFE post_prelimWorkMaps post_prelimFIRMs trtXpost_newmap"
		estpost summarize `variables' if lotFEsample==1, listwise
		esttab using "output_materials/Tables/Descriptive_stats_lotFE.tex", cells("sum mean(fmt(2)) sd min max") substitute(sum Count mean Mean sd Stdev min Min max Max) ///
			label nomtitle nonumber replace
			


*** Table A5 - With HFIAA

		eststo clear
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' `hfiaa'  if nbrhdFEsample==1 & blockFEsample==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' `hfiaa'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' `hfiaa' if lotFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Lot FEs -  
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' `hfiaa'   if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
		* Table
		esttab using "output_materials/Tables\DDmodels_heterog_dmgcontrol_hfiaa.tex", keep(1.flood_07#1.post_BW12 *floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy  1.postHFIAA#1*) ///
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy" 1.postHFIAA#1.flood_13_6 "HFIAA") ///
			mtitles("Neighborhood FE" "Block FE" "Block FE" "Lot FE") nonotes ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace	
			
		 

*** Table A6 - Linear Trends	by geography

		local cntrls "i.block i.yrqtr i.block#i.yrqtr" 

		eststo clear
		* Block & Lot FE (repeat columns 2 & 4 of primary table)
			eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
			eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)

		* Neighborhood linear trends by sale week
			eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols' i.neighborhoodINT##c.yrwk) vce(cluster census_tract_INT)
			eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols' i.neighborhoodINT##c.yrwk) vce(cluster census_tract_INT)

		*** Neighborhood by sale week FEs
			eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols' i.neighborhoodINT##c.yrwk) vce(cluster census_tract_INT)
			eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols' i.neighborhoodINT#i.yrwk) vce(cluster census_tract_INT)

		* Tax-block linear trends by sale week
			eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols' i.blockINT##c.yrwk) vce(cluster census_tract_INT)
			eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols' i.blockINT##c.yrwk) vce(cluster census_tract_INT)

		* Tax-lot linear trends by sale week
			eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols' i.lotINT##c.yrwk) vce(cluster census_tract_INT)
			eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols' i.lotINT##c.yrwk) vce(cluster census_tract_INT)

		* Table
		esttab using "output_materials/Tables/DDmodels_TaxLotLinearTrends_updated.tex", keep(1.flood_07#1.post_BW12 floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy) ///
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			mtitles("Block FE" "Lot FE" "Block FE"  "Lot FE" "Block FE"  "Lot FE" "Block FE"  "Lot FE" "Block FE"  "Lot FE" ) nonotes ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace

	
*** Table A7
		*** Generate an indicator for being "close to" the inundation area
		capture drop dist_to_flood
		gen dist_to_flood= (flood_sandy==0) * dist_2007_FP_edge/100
		generate nearby = inrange(dist_to_flood, .25, 1) /* .25 excludes properties very close to inundation to allow for measurement error in inundation boundary */
		gen logic = 0.flood_sandy#1.flood_07#1.nearby#1.post_Sandy /* unflooded, old FP, near inundation, after Sandy*/



		eststo clear	
		*** Baseline
		*** Block FEs - 13,490 in sample
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 i.flood_sandy `sandy2' alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)

		*** Lot FEs - 510,615 in sample
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)	


		*** Interact distance to flood with BW treatment			
		*** Block FEs - 13,490 in sample
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07 i.flood_07##i.post_BW12 1.flood_07#1.post_Sandy#c.dist_to_flood i.flood_sandy `sandy2' 1.logic 1.logic#c.dist_to_flood alltreated `newmaps' if nbrhdFEsample==1 & blockFEsample==1 , absorb(blockINT `timecontrols') vce(cluster census_tract_INT)

		*** Lot FEs - 510,615 in sample
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 1.flood_07#1.post_Sandy#c.dist_to_flood `sandy2' 1.logic 1.logic#c.dist_to_flood `newmaps' if lotFEsample==1 , absorb(lotINT `timecontrols') vce(cluster census_tract_INT)			

								
		*Output the Table
		esttab using "output_materials/Tables/DDmodels_BW_Dist_Fld_Table.tex", keep(1.flood_07#1.post_BW12 1.logic 1.logic#c.dist_to_flood floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth* trtXpost_newmapX*flood_sandy ) ///
			label numbers nonotes  collabels(none)  ///
			coeflabels(1.logic "Post-Sandy*in old FP*unflooded" 1.logic#c.dist_to_flood "Post-Sandy*in old FP*dist to fld" 1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			order(1.flood_07#1.post_BW12 1.flood_07#1.post_Sandy#c.dist_to_flood floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth* 1.logic 1.logic#c.dist_to_flood trtXpost_newmapX*flood_sandy) ///
			mgroups( "\underline{Main Estimates}"     "\underline{Hetero by Dist}"  "\underline{Full Interaction}" , pattern(1 0 1 0 1 0) span prefix(\multicolumn{2}{c}{) suffix(})) ///
			mtitles( "Block FE" "Lot FE" "Block FE" "Lot FE" "Block FE" "Lot FE" "Block FE" "Lot FE") ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace		
		


*** Table A8 - Restricting Sample Period

		*** Set Samples
				quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1 & inrange(saleyear,2012,2015), absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
						gen nbrhdFEsample_1215 = e(sample)
						quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1 & inrange(saleyear,2012,2015), absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
						gen blockFEsample_1215 = e(sample)
						quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1 & inrange(saleyear,2012,2015), absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
						gen lotFEsample_1215 = e(sample)
				
		*** Set samples
				quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1 & (post_Sandy!=1 | saleyear>2014), absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
						gen nbrhdFEsample_no1314 = e(sample)
						quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1 & (post_Sandy!=1 | saleyear>2014), absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
						gen blockFEsample_no1314 = e(sample)
						quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1 & (post_Sandy!=1 | saleyear>2014), absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
						gen lotFEsample_no1314 = e(sample)
								

		* Only 2012-2015
		eststo clear
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'  if nbrhdFEsample_1215==1 & blockFEsample_1215==1 , absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample_1215==1 & blockFEsample_1215==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)

			

		* omit post-Sandy Period
		eststo clear 
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'  if nbrhdFEsample_no1314==1 & blockFEsample_no1314==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample_no1314==1 & blockFEsample_no1314==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)


		* Table
		esttab using "output_materials/Tables\DDmodels_heterog_dmgcontrol_time_res_samples.tex", keep(1.flood_07#1.post_BW12 *floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy   ) /// 
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			mtitles("\shortstack{2012-2015\\ Neighborhood FE}"    "\shortstack{2012-2015\\ Block FE}" "\shortstack{Post-Sandy Omitted\\ Neighborhood FE}"    "\shortstack{Post-Sandy Omitted\\ Block FE}") nonotes ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace	


*** Table A9 - Future Floodplains
		* Calculating avg obs per cluster 
		preserve
		keep if nbrhdFEsample==1 & blockFEsample==1
		collapse (count) obs_within_cluster=logsaleprice_2010USD, by(census_tract_INT)
		display "Obs per cluster, columns 1-2"
		summarize obs_within_cluster
		restore
		preserve
		keep if lotFEsample==1
		collapse (count) obs_within_cluster=logsaleprice_2010USD, by(census_tract_INT)
		display "Obs per cluster, columns 3-4"
		summarize obs_within_cluster
		restore

		local future_2020_1 "i.in_2020_1##i.post_PlaNYC i.in_2020_1" 
		local future_2050_1 "i.in_2050_1##i.post_PlaNYC i.in_2050_1"


		* Adding flood damage to interacted model
		eststo clear
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' `future_2020_1' `future_2050_1' if nbrhdFEsample==1 & blockFEsample==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' `future_2020_1' `future_2050_1'  if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'  `future_2020_1' `future_2050_1' if lotFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Lot FEs -  
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps'  `future_2020_1' `future_2050_1'  if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
		* Table
		esttab using "output_materials/Tables/DDmodels_heterog_dmgcontrol_Fut_FPs.tex", keep(1.flood_07#1.post_BW12 *floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth   trtXpost_newmapX*flood_sandy 1.in_2020_1#1.post_PlaNYC 1.in_2050_1#1.post_PlaNYC) ///
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' 1.in_2020_1#1.post_PlaNYC "2020 1pct Floodplain" 1.in_2050_1#1.post_PlaNYC "2050 1pct Floodplain" trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			mtitles("Neighborhood FE" "Block FE" "Block FE" "Lot FE") nonotes ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace		
			


*** Table A10 - Hedonic Controls and Winsorized Prices	
		*** With Hedonic Controls
		eststo clear

		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'  `hed_cntrls'  if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)

		*** Lot FEs -  
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps'  `hed_cntrls'  if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)


		*** Winsorize, rather than trim outliers
			preserve
				winsor logsaleprice_2010USD, gen(win_logsaleprice_2010USD) p(0.01)
				replace logsaleprice_2010USD=win_logsaleprice_2010USD
				drop *FEsample
			
			* Quietly running models to generate sample flags with winsorized outcome 
				quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
				gen nbrhdFEsample = e(sample)
				quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
				gen blockFEsample = e(sample)
				quietly reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 i.flood_sandy#i.post_Sandy 1.trtXpost_newmap if inrange(saleprice, 100000, 6750000) & taxclassattimeofsale==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
				gen lotFEsample = e(sample)
				
				
				* Rerun saturated main models on winsorized outcome (and adjusted samples)
				*** Block FEs -  
				eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps'   if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
				*** Lot FEs -  
				eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps'  if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
				
			restore

			* Table
		esttab using "output_materials/Tables\DDmodels_heterog_dmgcontrol_hed_cntrls.tex", keep(1.flood_07#1.post_BW12 *floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy   ) /// 
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy") ///
			mtitles("Block FE" "Lot FE" "Block FE" "Lot FE") nonotes ///
			mgroups( "\underline{Additional Hedonic Controls}"     "\underline{Winsorized Sale Prices}"  , pattern(1 0 1 0 ) span prefix(\multicolumn{2}{c}{) suffix(})) ///
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace		
			

*** Tables A11 & A12 - Output of individual calculations entered by hand

		
**** Table A13 - Consider Flood protective infrastructure *	
		

		*** BigU Only	
		* Adding flood damage to interacted model
		eststo clear
		*** Mo, yr, neighborhood FEs
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' i.in_big_u i.in_big_u#i.post_big_u if nbrhdFEsample==1 & blockFEsample==1, absorb(neighborhoodINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' i.in_big_u i.in_big_u#i.post_big_u if nbrhdFEsample==1 & blockFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Block FEs -  
		eststo: reghdfe logsaleprice_2010USD `xsec' i.flood_07##i.post_BW12 `sandy2' `newmaps' i.in_big_u i.in_big_u#i.post_big_u if lotFEsample==1, absorb(blockINT `timecontrols') vce(cluster census_tract_INT)
		*** Lot FEs -  
		eststo: reghdfe logsaleprice_2010USD i.flood_07##i.post_BW12 `sandy2' `newmaps' i.in_big_u i.in_big_u#i.post_big_u if lotFEsample==1, absorb(lotINT `timecontrols') vce(cluster census_tract_INT)
		* Table
		esttab using "output_materials/Tables/DDmodels_w_Big_U.tex", keep(1.flood_07#1.post_BW12 floodsandyXpostsandyX*flood07 1.floodsandyXpostsandyX*flood07#c.sandy_inund_depth trtXpost_newmapX*flood_sandy  1.in_big_u#1.post_big_u ) ///   
			coeflabels(1.flood_07#1.post_BW12 "Biggert-Waters" `sandy2labels' trtXpost_newmapX1flood_sandy "Floodplain maps*Sandy" trtXpost_newmapX0flood_sandy "Floodplain maps*no Sandy" 1.in_big_u#1.post_big_u "Big U Protection" ) /// 
			mtitles("Neighborhood FE" "Block FE" "Block FE" "Lot FE") nonotes /// 
			se star(* 0.10 ** 0.05 *** 0.01) compress tex replace	
	

*******************************************************************************
*** Conduct Oster Unobservables Bounding Analysis in Section 7.1	
*******************************************************************************

		* Oster estimator can only use absorbed FE after areg, not xtreg
		* Her estimator doesn't use SEs, so not clustering to speed computation
		* R^2 in most saturated model is .84, and Oster's suggestion is 1.3*R~=1.09. That means default Rmax of 1 is probably appropriate here.
		* Oster-written command is "psacalc"
		* psacalc beta: estimates treatment effect assuming selection on observables and unobservables equal
		* psacalc delta: estimates relative selection on unobservables required to get zero treatment effect

		* Need to build local for psacalc
		* Allows us to focus on cross-sectional selection
		local psaoption ""
		forval i=2/763 {
			local psaoption "`psaoption' `i'.yrwk"
		}
		display "`psaoption'"

		* Tax lot FE
		areg logsaleprice_2010USD i.yrwk  1.flood_07#1.post_BW12 `sandy2' `newmaps' if lotFEsample==1, absorb(lotINT)
		*psacalc delta trtXpost_newmapX0flood_sandy, mcontrol(`psaoption')
		psacalc delta trtXpost_newmapX0flood_sandy

			


	
* Closing log
capture log close


				
